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Abstract. The running coupling and the Kugo-Ojima parameter of un¬ 
quenched lattice Landau gauge are simulated and compared with the contin¬ 
uum theory. Although the running coupling measured by the ghost and gluon 
dressing function is infrared suppressed, the running coupling has the maxi¬ 
mum of «o ~ 2 — 2.5 at around q = 0.5GeV irrespective of the fermion actions 
(Wilson fermions and Kogut-Susskind(KS) fermions). The Kugo-Ojima pa¬ 
rameter c which saturated to about 0.8 in quenched simulations becomes con¬ 
sistent with 1 in the MILC configurations produced with the use of the Asqtad 
action, after averaging the dependence on polarization directions caused by 
the asymmetry of the lattice. Presence of 1 + c\/q 2 correction factor in the 
running coupling depends on the lattice size and the sea quark mass. In the 
large lattice size and small sea quark mass, ci is confirmed of the order of a 
few GeV. The MILC configuration of a = 0.09fm suggests also the presence 
of dimension-4 condensates with a sign opposite to the dimension-2 conden¬ 
sates. The gluon propagator, the ghost propagator and the running coupling 
are compared with recent pQCD results including an anomalous dimension 
of fields up to the four-loop level. 


1 Introduction 

In 1978 Gribov showed that the Landau gauge fixing of the Yang-Mills theory 
does not define the unique gauge but there appear gauge equivalent copies and 
that the sufficient condition for the color confinement is the infrared vanishing of 
the gluon propagator [I]. A sufficient condition for the color confinement based 
on the Lagrangian field theory satisfying the Becchi-Rouet-Stora-Tyutin(BRST) 
symmetry was also proposed by Kugo and OjimafZj. Zwanziger showed that the 
uniqueness of the gauge field can be achieved by restricting the gauge field to 
the fundamental modular region and showed a horizon condition, which can be 
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checked by the lattice simulation [3j Investigations of the color-confinement 
criterion in Landau gauge as well as in the Curci-Ferrari gauge based on the 
Dyson-Schwinger approach are referred in ref. [5j. 

Lattice Landau gauge QCD simulation is a valuable tool for analyzing con¬ 
finement and chiral symmetry breaking from the first principle. In the simula¬ 
tion of quenched lattice Landau gauge of f5 =6, 6.4 and 6.45 with lattice volume 
24 ,32 4 ,48 4 and 56 4 , presence of infrared fixed point of «o ~ 2.5 was suggested 
and the tendency of increasing Kugo-Ojima parameter c as the continuum limit 
is approached was observedf6[[7]- The parameter c which is expected to be 1 for 
the proof of the confinement remained about 0.8. The lattice Landau gauge QCD 
suffers from gauge non-uniqueness problem i.e. Gribov copy problemjL, 8j and 
in large lattice we observed exceptional samples which possess axes along which 
the reflection positivity of the 1-d Fourier transform is violated and whose Kugo- 
Ojima parameters are close to 1. Most of the 1-d Fourier transforms of gauge 
configurations in quenched simulation violate rotational symmetry, but coupling 
to fermions recovers the symmetry. Besides this feature, a certain light meson 
propagator in quenched theory exhibits chiral loop artefacts|9’ [IQ], and thus 
the unquenched simulation results could be qualitatively different from those of 
quenched simulation. 

Recently, the MILC collaboration claimed that lattice QCD with three flavors 
agrees with variety of quantities with both light (u, d and s) and heavy (c or b) 
quarks with errors of 2-3%, whereas quenched QCD has errors as large as 15- 
20%[LI]. The Kogut-Susskind(KS) fermion contains artificial flavor degrees of 
freedom and MILC collaboration eliminated these degrees of freedom by taking 
the 4th root of the fermion determinant. This procedure can be justified when 
the flavor (taste) symmetry is preserved, but in the usual KS fermion approach it 
is violated in square order of the lattice spacing. The Asqtad action used by the 
MILC collaboration has an advantage that the taste violation would be reduced 
when the lattice spacing a is small [12] . We need to check this by measuring the 
running coupling and compare with the results of other fermions such as Wilson 
fermions, where no taste problem exists. 

We investigated gauge configurations of unquenched simulation in the data 
base i.e. JLQCD|13], GP-PACS|T4]. MILC [15] and Columbia University (CU)[l6l 
□2]. JLQCD and CP-PACS use Wilson fermions. The former is based on 0(a) 
improved Wilson action with a non-perturbatively defined clover coefficient csw ■ 
and the latter is based on the Iwasaki-improved gauge action with the tadpole 
improved clover coefficient csw- MILC and CU use KS fermions. MILC is based 
on the Asqtad improved action i.e. an extension of the Liischer-Weisz improved 
gauge action |I8, T9] and tadpole-improved fermion action. CU is based on the 
old standard Wilson action. 

In chiral perturbation theory, the length scale L = U 1 / 4 , the pion mass m^, 
infinite-volume pion decay constant F, quark condensate A and the effective 
cutoff A ~ 47 tF characterize the system. In order that particles fits well inside 
the box one requires the Compton wavelength of the pion l/m ff <C L. On the 
other hand, near the chiral limit extrapolation to 1 /m n ^$> L is required so that 
the collective Goldstone boson can be properly taken into account[22], In the case 
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of the Wilson fermion, there appears a problem to bridge the two regions due to 
appearance of a parity- and isospin-violating Aoki phase, and improvement by 
the twisted mass fermion etc. is proposed[23]. Whether KS fermion suffers from 
the same problem is discussed by several authors|24]. Thus it is important to 
clarify the infrared features of Wilson fermions and KS fermions on the lattice. 

The gauge configurations that we investigate are summarized in Tablc{TJ 


Table 1. j3,K sea and the sea quark mass m VWI (vector Ward identity) and the inverse lattice 
spacing 1/a, lattice size and lattice length(fm). Suffices c and / of MILC correspond to coarse 
lattice(a=0.12fm) and fine lattice(a=0.09fm). (3im P = 5/3 x /3. 
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Ksea 

VWI / VWI 

am ud / am s 

N f 

l/a(GeV) 

L s 

L t 

aL s (fm) 

JLQCD 

5.2 

0.1340 

0.134 

2 

2.221 

20 

48 

1.78 


5.2 

0.1355 

0.093 

2 

2.221 

20 

48 

1.78 

CP-PACS 

2.1 

0.1357 

0.087 

2 

1.834 

24 

48 

2.58 


2.1 

0.1382 

0.020 

2 

1.834 

24 

48 

2.58 

CU 

5.415 


0.025 

2 

1.140 

16 

32 

2.77 


5.7 


0.010 

2 

2.1 

16 

32 

1.50 

MILC c 

6.83 (flimp) 


0.040/0.050 

2+1 

1.64 

20 

64 

2.41 


6.76 (Pimp) 


0.007/0.050 

2+1 

1.64 

20 

64 

2.41 

MILC/ 

7.11 (Pimp) 


0.0124/0.031 

2+1 

2.19 

28 

96 

2.52 


7.09(f3imp) 


0.0062/0.031 

2+1 

2.19 

28 

96 

2.52 


In Sect.2 the unquenched lattice-simulation method is summarized and in 
Sect.3 numerical results of the gluon propagator, ghost propagator, Kugo-Ojima 
parameter and the running coupling are given. The conclusion and a discussion 
are presented in Sect.4. The pQCD formulae[39, 35] that are used in the fit of 
the lattice data are summarized in the Appendix. 


2 Unquenched Lattice simulation 


In the present lattice simulation, we adopt the log U type gauge field definition: 


Ux,fi — 


— pT,)! /it 
° ) X,!X 


-A 


X,fl- 


(1) 


The Landau gauge, dA 9 = 0 is specified as a stationary point of some opti¬ 
mizing functions Fjj(g) along gauge orbit [25] [3j where g denotes gauge transfor¬ 
mation, i.e., 5Fu(g) = 0 for any 5g. 

Here Fjj{g) for this options is [26, 6] 

Fu(s) = WII 2 = £*■' , (2) 

X,[l 


Under the infinitesimal gauge transformation g 1 5g = e, its variation reads for 
this defintion 


AFu(g) = -2{8A 9 \e) + (e| - 8D(U 9 ) |e) + • • • 
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where the covariant derivativative D^U) reads 



D^U^cj) = S(u x ^)d^ + [A x ^ $\ 

(3) 

Here 

= (j)(x + n) , &nd<j) = 

(4) 

and the definition of operation S(U X ^)B X41 is given as 



S{U Xjll )B x = T{A x ^)B x>tl 

(5) 

where A X4 <, 

x / 2 

= adjA m>ll = [Ac#, ■], and T(x) = th(r y 9) • 



The gauge fixing of the unquenched configuration is essentially the same as 
that of the quenched configuration. The convergence condition for the conjugate 
gradient method is less than 5% in the L 2 norm. We measure the gluon propaga¬ 
tor, the ghost propagator, the running coupling and the Kugo-Ojima parameter. 
Fermions will affect the running coupling through quark condensates, if they 
exist and will indirectly affect the Kugo-Ojima parameter. 

3 Numerical results 

In this section we show the lattice results and compare the data with the con¬ 
tinuum theory based on the effective charge method(27, 28]. In this method, the 
Green function G that depends on the scheme and the scale // is a solution of 
the renormalization group equation 

^(h, /i) = 7 (h)G(h, fj) (6) 

and the scheme and scale invariant G is expressed as 

G = G(h,fj,)/f(h) ( 7 ) 

where 

„. , f h dx 7 (x ). 

/(ft) = exp{j -± 1 ). ( 8 ) 

The function 7 ( 2 ;) is a function of the coupling constant and not a function of 
the scale.@In the MOM scheme, fi should be chosen such that the G(h,/j,) can 
be factorized into the scale-dependent part and the cut-off-dependent part. We 
have chosen y = 1.97GeV which is the inverse lattice spacing of the quenched 
(3 = 6.0 configuration and expressed propagators in terms of the effective coupling 
constant hm . The effective coupling constant is expanded by a parameter y which 
is a solution of 

I = A log(7M ~) - j a MM- (9) 

where A is the cut-off in the MOM scheme, (3q and (3\ are scheme inde¬ 
pendent constants of the /? function. A similar analysis using the principle of 
minimal sensitivity[29] (PMS) was performed in |30j , in which n is not fixed and 
the optimal parameter y is searched for each /i. The fit of the ghost propagator 
in the PMS was not successful, and we adopted the MOM scheme. Details of 
the effective charge method are summarized in the Appendix. 
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In the log— U definition we express the gauge field defined at the midpoint of 
the link as 

U x ,n = exp A X4l = exp[y A* [x + ^)] (10) 

A a ^ X b 5 ab 

where A a is the Gell-Mann’s A normalized as tr-= —. 

2 2 2 

This definition is consistent with the continuum theory [21]. but differ from 
the usual lattice convention [20] 

[/;(*) = ex p[AM“]. (11) 

where the lattice spacing is taken as 1 and the coupling constant is included in 

The gluon propagator Da(q 2 ) in the Landau gauge is defined as 


Aa 

/w 


D A,fii/(Q ) — 


m -1 


tr(Ajq)Al(q)) 




( 12 ) 


where N c is the number of color and the Fourier transform of the gauge field is 

= (13) 


The summation over p gives 

N 2 _ x t r (i M (g)^(g)) = (d- 1 )D A {q 2 ) (14) 

c a 

The midpoint definition is crucial in the definition of the momentum of the 
Liischer-Weisz’s improved action. In the practical calculation, however, one is 
allowed to perform the simple numerical Fourier transform 

Mi) = 

= <?*'*''Mi) < 15 ) 

and evaluate 

N 2 _ 1 M^(<?)^(<?)) = (d - i )D A (q 2 ) (16) 

c a 

which yields trivially the same Da{q 2 ) as that obtained by the A^(q). 

In passing, we remark that the Ansatz of the gauge transformation adopted 

by [32] 


% 0*0 = e eX ^ x %(x) 


( 17 ) 
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where in the classical continuum limit X p (x) = —a 3 '^^D 1/ F IJjI/ (x) yields trans- 

V 

formation 

Afj, + ^a 2 ^(d v F pv (x) + [A v , F pv (x)}) + o(a 4 ) (18) 


In our Ansatz 

Upix) = glu p (x)g x+p (19) 

where g\ = exp[|a 2 ^ F^ v {x)\ and g x+p = exp[^a 2 ^ F pu (x + //)], yields the 

~ V V 

same expressioin of o(a 2 ) in the Landau gauge. In the derivation of the Liischer- 
Weisz improved action, the definition U p r(x) = exp aA^it7(x) was used[18], 
where R is an irreducible representation of SU(N). 

In the data analysis of we usually adopt q diagonal in the momen¬ 

tum lattice, which is called cylinder cut. In the case of different lattice spacial 
length N s and time length Nt, we define cylinder cut as qa around the diagonal 
[qi,Q2,Q3, Qa] = [q, q, q, (Nt • q/N s )] where (Nt • q/N s ) is the closest integer to the 
quotient of Nt ■ q/N s . 

When the improved action is adopted, the lattice gluon momentum is defined 
asfHJS^ q = \Zq 2 + ^ 4 /3/a where 


3 


4<y; sm 2 (qiTr/N s ) + sin 2 (q ; 47r/IVi)) 

2—1 

(20) 

3 

4<y; sin 4 (qiir / N s ) + sin 4 (q^/N t )) 

2=1 

(21) 


/ q 

The correction factor W 1 H- tt is about 1.15 in the highest momentum point of 

V 3 r 

cylinder cut. 

In Fig. |U the logarithm of the gluon dressing function as a function of 
logio q(GeV) of the Wilson fermion K sea = 0.1357 and K sea = 0.1382 are shown. 
The corresponding data of the KS fermion f3i mp = 7.11 and 7.09 are shown in 
Fig. [2] In the infrared region the dressing function of light Wilson fermion mass 
(. K sea = 0.1382) is enhanced as compared to that of the heavy Wilson fermion 
mass ( K sea = 0.1357), while that of heavy KS fermion mass (f3i mp = 7.11) and 
that of light KS fermion mass (Pimp = 7.09) are almost identical. 

In Figs. [3] and [4] the gluon propagator and the gluon dressing function of 
MILC fine lattice (MILC/, Pimp = 7.09) in cylinder cut are shown. Data of 
MILC coarse lattice (MILC c ) are consistent with those of [33]. In [33], the gluon 
propagator is normalized as 1 /q 2 at q = 4GeV and the data are about factor 2 
smaller than ours. 

The pQCD formulae of the gluon dressing function are given in the Appendix. 
We observe that the gluon propagator of pQCD in quenched approximation is 
larger than that of the unquenched one and the result of the four-loop calculation 
is larger than that of the three-loop calculation. (Fig. [5]) The corresponding data 
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Log_10[q(GeV)] 


Figure 1 . The log of the gluon dressing function of CP-PACS K sea = 0.1357(diamonds)(50 
samples) and 0.1382(triangles)(50samples). 



-1 -0.5 0 0.5 

Log_10 [q(GeV)] 


Figure 2. The log of the gluon dressing function of MILC /3im P = 7.09(stars)(50 samples) and 
7.11 (triangles) (50samples). 


of gluon dressing function are shown in Fig. [ 6 l In these plots we used = 
0.237GeV, y = 1.97GeV, y = 0.0222703 and A = 17.85 to fit Z(9.5GeV)=1.3107 
[Mj obtained in the quenched Landau gauge simulation. The definition of the 
effective coupling strength y and the strength A are given in the Appendix. 

The difference of Z(q) between lattice data and pQCD results would yield 
information on the gluon condensates [40j. We performed x 2 fit of the difference of 
the gluon dressing function of MILCj Pi mp = 7.09 and pQCD four-loop Nf = 3 
result. We parametrize 


Ziattiq 2 ,/x 2 ) = Z pQCD (q 2 ,y 2 )( 1 + ^ )+d (22) 

where y = 8.77GeV and the data of q > 3GeV region fit by searching ci and d. 
Fig{T] is the fitting result using the pQCD result for y = 8.78GeV, y = 0.0148488 
which gives c\ = 7.39GeV 2 , d = —0.024, x 2 /d.o./ = 1.10 . The fit using y = 
1.97GeV gives c.\ = 7.04GeV 2 , d = —0.017 and x 2 /d.o./ = 1.14. 
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0 2 4 6 8 

q(GeV) 


Figure 3. The gluon propagator of MILC/ /3j mp = 7.09. 



q(GeV) 


Figure 4. Same as Fig. [3] but the gluon dressing function. 


3.2 The ghost propagator 

The ghost propagator is defined by the Fourier transform(FT) of the expectation 
value of the inverse Faddeev-Popov(FP) operator j\4 

FT[Dq(x, y)\ = FTitriA^UMiU})- 1 }^), 

= 5 ab D G (q 2 ). (23) 

where antihermitian A a is normalized as tryl“'A & = 5 ab . 

The ghost dressing function G(q 2 ) is defined as 

D G (q 2 ) = (24) 

r 

The pQCD formula of the ghost dressing function is given in the Appendix. 
We observe that the ghost propagator of pQCD in quenched approximation is 
larger than that of the unquenched one, and the result of the four-loop calculation 
is smaller than that of the three-loop calculation. (Fig. [8|) The corresponding data 
of the ghost dressing function are shown in Fig. [9j 

In these plots we used = 0.237GeV, p = 1.97GeV, y = 0.0158465 and 

A g = 1- 
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q(GeV) 

Figure 5. The gluon propagator of pQCD 3-loop Nf = 0(short-dashed),Af/ = 3(dash-dotted). 
4-loop Nf = 0(solid),Al/ = 3(dashed). 



q(GeV) 


Figure 6. Same as Fig. [5] but the gluon dressing function. 


The ghost propagator of CP-PACS K sea = 0.1382 is shown in Figs. [[Ill The 
solid curve in Fig. [TO] is the pQCD result in the three-loop Nf = 00 , using 
the scale parameter y = 1.97GeV and Xq = 3.22 that fit the quenched data. 
The dashed line is the pQCD result in the four-loop Nf = 2 using the same 
A g and y. In the region 2 < q < 7GeV the pQCD result of Nf = 2 is about 
10% smaller than that of Nf = 0 and parameters A q = 3.22, y = 0.0158465 give 
qualitatively good agreement. In order to get better agreement of CP-PACS, we 
change y = 0.024610 as the solution of Nf = 2 and perform y 2 fit of Xq using 
the data of ghost dressing function in q > 2.8GeV region. The fit with A q = 3.01 
is shown in Fig 1 111 

The logarithm of the ghost dressing function as a function of log 10 q(GeV) of 
the Wilson fermion (CP-PACS) and the KS fermion (MILC) are shown in Figsll2l 
and[[3l respectively. In the case of CP-PACS, the ghost dressing functions are 
almost independent of the K sea and the exponent olq in the region q > 0.4GeV 
is 0.22(5), and in the case of the MILC fine lattices it is 0.24(5) in Pi mp = 7.09 
and 0.23(5) in Pi mp = 7.11. The exponent of the quenched, the unquenched 
Wilson and the KS fermions are almost the same. In view of the discontinuity 
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0 2 4 6 8 

q(GeV) 


Figure 7. The difference AZ(q 2 ) = Zi at t(q 2 , y 2 ) — Z pert (q 2 , y 2 ) as a function of q( GeV). The 
dashed line is a fit by condensates of dimesion-2. 



q(GeV) 


Figure 8. The ghost propagator of pQCD 3-loop Nf = O(short-dashed), Nf = 3(dash-dotted). 
4-loop Nf = 0(solid),IV/ = 3(dashed). 


in the slope of the ghost propagator in the infrared region, we ignore the lowest 
few points in the evaluation of the running coupling. In the asymptotic region 
(log 10 q(GeV) ~ 1) the ghost dressing function of MILC c , MILC/ and CP-PACS 
converge as shown in Figsll2l andfT3l 

We performed x 2 fit of Xq for ghost dressing function of MILCj. Using y = 
1.97GeV and Nf = 3, we obtained y = 0.026775 but with this parameter the 
fit was not better than that using y = 0.024610 i.e. the solution of Nf = 2. We 
obtained Xq = 3.258, y = 0.024610 for (3i mp = 7.09 data and A g = 3.237, y = 
0.024610 for /3i mp = 7.11 data. Since the sea-quark mass rn s is large, the better 
fit using pQCD with Nf = 2 is not unexpected. It is remarkable that Xq of 
MILCj are close to that of quenched simulations, i.e. A g = 3.22. 
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Figure 9. Same as Fig. [8] but the ghost dressing function. 



q(GeV) 

Figure 10. The ghost propagator as a function of q(GeV) of CP-PACS K sea = 0.1382(50 
samples). Solid curve is the 3-loop Nf = 0 pQCD result and dashed curve is the 4-loop Nf = 2 
pQCD result (Ac = 3.22, y = 0.0158465). 


3.3 The Kugo-Ojima parameter 

The Kugo-Ojima parameter is defined by the two point function of the covariant 
derivative of the ghost and the commutator of the antighost and gauge field 



(25) 

The confinement criterion is that the parameter c defined as u ab (0) = —5 ab c 
becomes 1. The parameter is related to the renormalization factors which are 
defined in the next section as m 
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4 
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3 

__ 2 .5 

5 2 

U 1.5 
1 

0.5 


Figure 11. The ghost dressing function as a function of q( GeV) of CP-PACS K sea = 0.1382(50 
samples). Dashed line is the 4-loop Nf = 2 pQCD result (Ag = 3.01, y = 0.0246100). 


T 

I 



2 4 6 8 

q(GeV) 



Log_10[q(GeV)] 

Figure 12. The log 10 G(q) as a function of log 10 q( GeV) of CP-PACS K sea = 
0.1357(diamonds)(50 samples) and K 3£a = 0.1382(triangles)(50samples). 


If the finiteness of Z\ is proved, divergence of Z3 is a sufficient condition. If Z3 
vanishes in the infrared, Z\ should have higher order 0. 

From the investigation of the Gribov problem, Zwanziger proposed the hori¬ 
zon function [(2 |] 


v 

x,y 




1 

-dD 


(■ ~D V )A 



= G, v {q)5 ah 


(27) 


where e = (E tr(A a 'S(U x ^)A a )\ /{{Nl - l)V}, and N c = 3 for SU(3). The 
\ / 

horizon condition reads lim. G fJjfl (q)—e = 0, and the left-hand side of the condition 

is (-) + (d — l)c — e = (d — 1 )h where h = c — 1 and dimension d = 4, and 
\dJ d 
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Log_10[q(GeV)] 

Figure 13. The log 10 G(q) as a function of log 10 q(GeV) of MILC/ Pimp = 7.09(stars)(6 sam¬ 
ples) and MILC c pimp = 6.76(triangles)(50 samples). 


it follows that h = 0 implies the validity of the horizon condition, and thus the 
horizon condition coincides with Kugo-Ojima criterion provided the covariant 
derivative approaches the naive continuum limit, i.e., e/d = 1. 

Therefore the infrared behavior of the gluon propagator and the ghost propa¬ 
gator, as well as direct measurement of the Kugo-Ojima parameter are important 
for understanding the color-confinement mechanism. 

The Kugo-Ojima parameter c of the unquenched simulation depends on the 
direction of the polarization due to asymmetry of the lattice. When the polar¬ 
ization is in the spacial directions, c is consistent with 1 in most unquenched 
simulations. (See Tabic lT3l) 


Table 2. The Kugo-Ojima parameter for the polarization along the spacial directions 
c x and that along the time direction ct and the average c, trace divided by the di¬ 
mension e/d, horizon function deviation h of the unquenched Wilson fermion(JLQCD, 
CP-PACS), and KS fermion (MILC c ,CU,MILC/). The log U definition of the gauge 
field is adopted. 



K aea or p 

C-x 

Ct 

c 

e/d 

h 

JLQCD 

K sea =0.1340 

0.89(9) 

0.72(4) 

0.85(11) 

0.9296(2) 

-0.08(11) 


K sea =0.1355 

1.01(22) 

0.67(5) 

0.92(24) 

0.9340(1) 

-0.01(24) 

CP-PACS 

K sea =0.1357 

0.86(6) 

0.76(4) 

0.84(7) 

0.9388(1) 

-0.10(6) 


K sea =0.1382 

0.89(9) 

0.72(4) 

0.85(11) 

0.9409(1) 

-0.05(9) 

CU 

p =5.415 

0.84(7) 

0.74(4) 

0.81(8) 

0.9242(3) 

-0.11(8) 


P =5.7 

0.95(26) 

0.58(6) 

0.86(28) 

0.9414(2) 

-0.08(28) 

MILC c 

P =6.76 

1.04(11) 

0.74(3) 

0.97(16) 

0.9325(1) 

0.03(16) 


P =6.83 

0.99(14) 

0.75(3) 

0.93(16) 

0.9339(1) 

-0.00(16) 

MILC/ 

P =7.09 

1.06(13) 

0.76(3) 

0.99(17) 

0.9409(1) 

0.04(17) 


p =7.11 

1.05(13) 

0.76(3) 

0.98(17) 

0.9412(1) 

0.04(17) 
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3.4 The running coupling 

The calculation of the running coupling in the MOM scheme using the gluon 
dressing function and the ghost dressing function is discussed in [47]. We 
parametrize the gluon dressing function as 

Z R (q 2 ,p 2 ) = Z^{/3,p 2 )Z(q 2 ) (28) 

and the ghost dressing function as 

G R {q 2 ,p 2 ) = Z~\p 2 )G{q 2 ) (29) 

with the renormalization conditions 

Z R (p 2 ,p 2 ) = G R (p 2 ,p 2 ) = 1 (30) 

where p is an accessible scale of the lattice simulation [461. We choose p ~ 6GeV 
. Infrared properties of Z{q 2 ) and G(q 2 ) are defined as Z(q 2 ) oc ( qa)~ 2oiA and 
G(q 2 ) oc (qa)~ 2a °. 

We define the vertex renormalization factor Z\ as 

aR(p 2 )Z R (q 2 , p 2 )G R (q 2 , p 2 ) 2 

= esgph x pi) 

{(3, p) 

where the subscript R means ’’renormalized”, and 

®R{q 2 ) = OL R {p 2 )Z R (q 2 , p 2 )GR(q 2 , p 2 ) 2 (32) 

The MILC collaboration adopted the tadpole improvement in the generation 
of the gauge configuration. In [46] the gluon propagator is calculated as 

a 0 (A uv )a 2 D A (x,y) = D A {x,y) /ul P (33) 

where Da(x, y ) corresponds to the propagator in terms of the link matrices, and 
the ghost propagator as 

a 0 (A uv )a 2 D G (x, y) = D G (x, y)u 0 ,p (34) 

where uo,p = (P) 1 ^ 4 . 

The plaquette value {P) of MILC c is smaller than that of MILC f by a few 
percent and the correction in the FigHHby this renormalization is negligible. In 
the running coupling, the tadpole factors cancel out[46|. We modify the notation 
of eq. ([3TT) and measure running coupling a s (q ) as[38] 

= ~ a,(A uv )r 2im+2a °\ (35) 

47t Zf 

where q = y/q 2 + q 4 (S /3)/ a. (4 = 0 for ordinary action and 4=1 for the 
improved action.) 
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In the quenched simulation the Orsay group fitted the lattice data by 
&s,iatt(p) = cts,pert(^) + ci//x 2 with c\ ~ 0.65GeV 2 [35]. In the analysis of un¬ 
quenched Wilson fermion data[36l 157] . they fitted the lattice data in the form 


&s,latt (p) 


— 0^s,pert (/t) (1 H“ 



(36) 


with ci ~ 2.8(2)GeV 2 . 

In the quenched simulation we confirmed the correction term [7] and we stud¬ 
ied whether the same correction appears in the unquenched configuration. We 
define the scale of the running coupling by fitting Z\ such that the running 
coupling at q ~ 6GeV agrees with the pQCD results of Nj = 2(JLQCD,CP- 
PACS,CU) or Nf = 3(MILC). In the case of the improved action, there is an 
ambiguity due to the mismatch of the momentum of the ghost dressing function 
q/a and that of the gluon dressing function q. 


Table 3. The \/Z\ factor of the unquenched SU(3). 


config. 

heavy 

light 

N f 

comments 

JLQCD 

0.90(7) 

0.97(7) 

2 

Ksea = 0.1340,0.1355 

CP-PACS 

1.07(8) 

1.21(10) 

2 

Ksea = 0.1357,0.1382 

CU 

1.13(10) 

1.19(8) 

2 

p = 5.415,5.7 

MILCc 

1.49(11) 

1.43(10) 

3 

Pimp = 6.83, 6.76 

MILC f 

1.37(9) 

1.41(12) 

3 

Pimp = 7.11,7.09 


The magnitude of the running coupling a s (q) in the infrared is roughly pro¬ 
portional to the 1 / Z'f factor. 

The mass of the sea quark is relatively heavy in JLQCD. There is a deviation 
from the pQCD in the region q < 3GeV, but above 3GeV the deviation is within 
statistical errors. The CP-PACS configuration has a lower sea quark mass and 
an Iwasaki-improved action is used for the gauge action. The running coupling 
of CP-PACS is shown in Fig. [TT] The absolute value increases as the mass of 
sea quark decreases and the data of lightest quark mass suggest a maximum of 
a s ~ 2 — 2.5. 

We measured the running coupling of the KS fermion in small /3 (strong 
coupling region) using the gauge configurations of CU[T6l [17] , and the large (3 
using those of MILC. In contrast to the Wilson fermion, the absolute value of the 
running coupling is close to that of JLQCD K sea = 0.1355 and does not depend 
on the mass of the sea quark. 

The MILC collaboration improved the flavor symmetry violation in the KS 
fermion by choosing an appropriate improved fermion action which is called 
Asqtad action. The running coupling of MILC c and MILC/ are shown in Figs. 
m anddHl respectively. 

We observe that the absolute value is consistent with that of the CP-PACS 
of the smallest sea quark mass. The correction of c\/q 2 with c\ of the order of 
2.8GeV 2 observed by the Orsay group in the Wilson fermion exists also in the 
KS fermion of Asqtad action. 

In Figs. [T7] and [18] the difference of running coupling of lattice data and 
pQCD data, Aa s = a S) i att — ot SyPer t as a function of q( GeV) for JLQCD K sea = 
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Figure 14. The running coupling a s (q) as a function of log 10 g(GeV) of CP-PACS (5 = 2.1, 
24 3 x 48 lattice, K sea = 0.1357 (diamonds) and that of 0.1382 (triangles), (25 samples each). 
The long dashed line is the Dyson-Schwinger fit with ao = 2.5, the dash-dotted line is the 
pQCD result and the short dashed line is the pQCDx(l + c/q 2 ) correction. 

3 

2.5 
2 

— 1.5 

CO 

a 

l 

0.5 

-0.4-0.2 0 0.20.40.6 0.8 1 

Log_10 [q(GeV)] 

Figure 15. Same as Fig. 1141 but the data of MILC c (a = 0.12fm) Pimp = 6.76(triangles) and 
6.83(diamonds), (50 samles each). 



0.1355, CP-PACS K S ea = 0.1382 are plotted, respectively. The fitted curves are 
a s (q)ci/q 2 + cL where the parameters c\ and d are obtained by the x-square fit 
using data above ~ 3GeV (solid) and data above ~ lGeV (dashed), respectively. 
In these configurations and in MILC c /3i mp = 6.83, the parameter c\ obtained by 
the two fits agrees within the statistical errors. The parameter ci ~ 2.4(2)GeV 2 
of CP-PACS is 15% smaller than the Orsay fit of K sea = 0.15 Wilson fermion 
data. The parameter c\ ~ 1.9(3)GeV 2 of MILC c (3i mp = 6.83 is about 2/3 of 
the Orsay fit. The a s of JLQCD have c\ ~ 1.15(4)GeV 2 but above 3GeV it is 
consistent with pQCD. 

In the case of MILC c (3 = 6.76 and MILCj, Aa s cannot be fitted by the factor 
-i, therefore we fitted above lGeV region by the factor H—r + "T an d the 

qz ■ ^4 go 

overall shift d. In the theory of operator product expansion, is proportional 
to the gluon condensates (o 2 F 2 t .)[40l [49J and/or quark condensates ( mqq ), and 
C6 is proportional to condensates of dimension 6 like (g 3 f^F^F^F^) . We 
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Figure 16. Same as Fig. [14] but the data of MILC/ (a = 0.09fm) pi mp = 7.09(stars) and 
7.11 (diamonds) (15samples each). 


performed x 2 fit of c t (i=l,2,3) by either fixing C3 = 0 or without fixing. As 
shown in Table l3~4l we find c± ~ 4.2(l)GeV 2 , C2 ~ —2.3(2)GeV 4 , C3 = 0 or 
ci = 6.6(l)GeV 2 , C2 = —13(2)GeV 4 , C3 = 8(2)GeV 6 as an average of f3i mp =6.76, 
7.09 and 7.11 data. The fit of Aa s of MILC c j3i mp = 6.76 data is shown in Fig. 
[T9l and that of MILCy is shown in Fig. 1201 



Figure 17. Aa s of JLQCD K = 0.1355 as a function of q(GeV). The fitted lines have ci = 
1.12GeV 2 and ci = 1.19GeV 2 , respectively. 



q(GeV) 

Figure 18. Aa s of CP-PACS I\ = 0.1382 as a function of g(GeV). The fitted lines have 
ci = 2.57GeV 2 and ci = 2.38GeV 2 , respectively. 

Although statistics is not large, running coupling of CP-PACS and MILC 
show a deviation from pQCD at q ~ 3GeV region and the deviation suggests a 
presence of A 2 condensates, gluon condensates and/or quark condensates. The 
origin of the enhancement of the running coupling of MILCj in 3-5GeV region 
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12345678 

q(GeV) 


Figure 19. Aa s of MILC c /3 = 6.76 as a function of q( GeV). The fitted lines have ci = 
6.50GeV 2 ,C2 = —ll.TOGeV 4 , C 3 = 7.47GeV 6 and d=-0.29 (solid) and ci = 4.18GeV 2 ,C2 = 
—2.45GeV 4 and d = —0.014 (dashed). 



Figure 20. Aa 3 of MILC/ /3i mp =7.11(diamonds) and 7.09(stars) as a function of q(GeV). The 
solid line is the fit of pimp = 7.11(ci = 4.27GeV 2 , C 2 = —2.28GeV 4 ) and dashed line is the fit 
of p^ = 7.09 (ci = 4.29GeV 2 , c 2 = -2.30GeV 4 ). 


is the enhancement of the ghost propagator in qa > 1 region (Fig ll3l) . Although 
X 2 becomes smaller for C3 p 0, we need further study for verifying presence of 
condensates of dimension 6. 


Table 4. The yp fit of coefficients d in of the 
MILC/ Pimp =7.11 and 7.09. C 3 fixed to 0 and unfixed 
cases, d.o.f means the degrees of freedom. 


Pimp 

Cl 

C2 

C3 

d 

X 2 /d.o.f 

6.76 

4.18 

-2.45 

0 fix 

-0.0143 

0.95 


6.50 

-11.7 

7.47 

-0.0289 

0.44 

7.11 

4.27 

-2.28 

0 fix 

-0.0281 

2.6 


6.58 

-13.05 

8.47 

-0.0395 

1.7 

7.09 

4.29 

-2.30 

0 fix 

-0.0344 

3.6 


6.58 

-13.84 

9.27 

-0.0449 

2.6 


In ref. [32j, the running coupling of the MILC configuration was measured 
by using a perturbative expansion for the plaquette and upsilon spectroscopy 
to set the scale. They found 07,(8.2GeV)~ 0.214 at Nf = 2 + 1, m^/mp ~ 0.4. 
When we fix the scale of MILC/ /3i m p = 7-09 by the Nf = 3 pQCD result at 
q =6.84GeV, a s (6.84GeV)=0.219, our data o s (8.2GeV)=0.190(2) is smaller than 
[42] by about 10%. A possible origin of the difference is that the pQCD results 
correspond to those of the chiral limit and does not fix the proper scale of the 
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lattice. A more rigorous scale fixing would be achieved by fitting lattice data in 
high energy region in perturbation series of the running coupling, fixing the a s 
in e.g. V-scheme and converting to the MS scheme {421 fl3l 144 1 . 

The V-scheme would not be applicable to the infrared region due to the 
sensitivity to the test-charge wave function. Brodsky et al.[45| applied physical 
a T scheme in hypothetical r lepton decay and, by using the /3 T three-loop results, 
observed freezing of the running coupling to the infrared fixed point a T (0) ~ 2. 


4 Conclusion and discussion 


We measured the running coupling and the Kugo-Ojima parameter from un¬ 
quenched QCD gauge configurations of Wilson fermion and KS fermion. We ob¬ 
served a sign of different behaviors of the Wilson fermions and the KS fermions 
as the system approaches to the continuum limit and the chiral limit. In the 
case of the Wilson fermion, the running coupling increases as the mass of the sea 
quarks decreases and approaches to the chiral limit, while in the case of the KS 
fermion dependence on the mass of the sea quark is very weak but it depends on 
the lattice spacing a. The a dependence of the KS fermion was expected to be 
due to the presence of violation of the taste symmetry which is of order a 2 and 
the Asqtad action improved this deficiency. Despite these differences, the run¬ 
ning coupling of the Wilson fermion of the smallest quark mass i.e. CP-PACS 
(/3 = 6.85) configuration and that of the KS fermion of the MILC/ and MILC c 
{Pimp = 6.76) are consistent. Milder finite size effects in Wilson fermions than 
in KS fermions, due to the spread of the KS fermion over a hypercube in the 
spinor-flavor interpretation is observed also in [48 ]. 

The ci / q 2 correction of the running coupling of ci of the order of a few GeV 
is confirmed in the CP-PACS data but not in JLQCD. Orsay group analyzed 
the data of Wilson fermion of K sea = 0.15, while we analyzed JLQCD Wilson 
fermion of K sea = 0.1340 and find the correction is much smaller. These results 
indicate that the term c\/q 2 appears as the system approaches the chiral limit. 
The MILC/ data suggests that near the chiral limit there are c\/q 2 and C 2 /q 4 
terms which have different signs. The different sign of the subsequent terms 
causes worriying whether this expansion converges. We need to increase the 
statistics for obtaining a definite magnitude of the condensates. 

On the physical meanings of the c\/q 2 term, there are several discus¬ 
sions. The Orsay group interprets this term as an indication of the A 2 
condensatesf35j 00) 09] • The operator A 2 is a dimension-2 operator allowed to 
have a vacuum expectation value, and it appears in the operator product expan¬ 
sion of the running coupling and in the gluon dressing function. Although it is 
not gauge invariant, the Landau gauge condition = 0 is compatible with 

stationarity of (A 2 )[50]. In the context of the maximal abelian gauge, the on- 
shell BRST invariant mixed gluon-ghost condensate of dimension-2 is discussed 
as gauge invariant observable [5T, 52]. At tree level the parameter ci in the gluon 
dressing function and ci in the running coupling in the triple gluon vertex are 
related to the (A 2 ) as[40] 


ci 


3 g z 


(A 2 ) 


prop 


4 (iV 2 - 1) 


(37) 
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and 


ci = 9 g' 


3 alpha 
- 1 ) 


(38) 


In the ghost anti-ghost gluon vertex, the multiplicity of ci is reduced by a factor 
3. Thus the fit of running coupling and the gluon dressing function with [i = 

^ 2 ) 1/2 

1.97GeV, C3 = 0 ansatz yields - al v^ a _ q 73 an( ^ the corresponding case 

{A 2 )prop 

with C3 / 0 ansatz yields 0.97. Orsay group obtained the ratio in the quenched 
simulation as 1.21(10] in the 3-loop calculation. 

In the restriction of the gauge field of the Landau gauge in the fundamental 
modular region, Zwanziger[53] defined the horizon function as ( H) = V(N 2 — 
1 ){dh + e), where V is the lattice volume and d,e and h are defined following 
Eg. (1271) . He proposed a simulation with a Boltzmann weight of e~ lH , where 7 1 / 2 
is a parameter of dimension-2. The dimensional parameter in action breaks the 
dilatation invariance and it has a link to the global properties of the fundamental 
modular region. This conjecture was recently discussed in ref. [54] including a 
possible condensation of A 2 . In this theory Zwanziger’s 7 is affected by the 
presence of A 2 condensates. Since our simulation is done without the Boltzmann 
weight, we cannot measure the parameter 7 directly. We observed, however, the 
horizon function h is negative in quenched simulations and consistent with 0 in 
unquenched large lattices simulations. The running coupling of the quenched as 
well as unquenched simulation shows that {A 2 ) is positive. The running coupling 
of MILC/ suggests that C2 where the gluon condensates ( g 2 F 2 u ) and/or the 
quark condensates ( mqq ) contribute is negative. In an analysis of QCD gap 
equation [54], a solution with negative ( g 2 F 2 v ), positive {A 2 ) was found. Since 
the gluon condensates ( g 2 F 2 u ) and the vacuum energy E vac have opposite signs, 
it implies that the sign of vacuum energy is positive in contradiction to the result 
of a two-loop analytical calculation [55]. 

Since C2 term appears only in the running coupling of MILC/ and MILC c of 
light sea quark mass, it would be natural to interpret that the C2 term comes 
mainly from quark condensates. Although the sea-quark mass of MILC/ m u + 
rrid = 2 X 0.068GeV is too heavy to discuss the chiral limit, (qq) has the correct 
sign as the Gell-Mann, Oakes, Renner relation [56l 57j 


m n 2 fl ~ -(m u + m d )(qq) + 0(m 2 ud ) (39) 


requires. An analysis of quark propagator [58] in Landau gauge also suggests that 
it is negative. 

Momentum dependence of the running coupling in the infrared region is 
chracterized by the infrared exponent of the gluon dressing function and the 
ghost dressing function. The infrared exponent ac of the ghost dressing function 
at 0.4GeV region is about half of k used in the Dyson-Schwinger approach at the 
infrared limit. We observe 2 «g + Q£) ~ 0 which supports the presence of infrared 
fixed point[63]. In the asymmetric lattice we observed that ac near the lowest 
momentum along the long lattice axis is suppressed. As an analysis of Dyson- 
Schwinger equation suggests, this suppression coud be due to the compactness 
of the lattice[63J - We suspect, however, there are effects due to the fluctuation 
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of the ghost propagator[55]. and there are an elaborated structure of the fixed 
point which is veiled by Gribov copies. 

Zwanzigerjl] argued that stochastic gauge fixing would render configurations 
to the common boundary of the fundamental modular region and the Gribov 
region and the Gribov copy effects can be evaded.@The argument is based on 
the assumption that the renormalization group flow of the ghost propagator 
follows perturbative renormalization-group flow equation. Our simulations of the 
ghost propagator [65] do not confirm the simple renormalization-group flow. The 
scenario of suppression of the infrared modes of the gauge held due to vanishing 
gluon propagator is not confirmed neither and yet to be further investigated. 

We admit that we could not restrict our gauge fixed configurations of un¬ 
quenched simulations in the fundamental modular region. In the case of SU(2), 
we performed parallel tempering(PT) gauge fixing and observed that the ghost 
propagator of PT gauge fixed samples is less singular than that of the first 
copy[6, 59]. In SU(3) of large lattice volume, we found exceptional copies whose 
A 2 norm is larger than the average and the ghost propagator is more singular 
than the average. These data suggest that infrared features of Gribov copies are 
complicated. 

Freezing of the running coupling in the infrared is assumed in a model of dy¬ 
namical chiral symmetry breaking [60l [61]. A Dyson-Schwinger approach predicts 
that the infrared fixed point corresponding to k = 0.5 is about 2.5ED and the 
lattice results would not be inconsistent with this model, if the lattice artefacts 
could be properly removed. 

The Kugo-Ojima parameter of MILC configuration is consistent with 1 in 
the average of polarization, but the value for polarization in t direction is small, 
since the lattice length transverse to t direction is short in asymmetric lattices. 
The slope of the ghost propagator ac also depends on the length of the axis. 
Differences between two directions provide a warning on lattice artefacts in the 
infrared region, however the qualitative difference of the quenched simulation 
(c ~ 0.8) and the unquenched (c ~ 1) would be related to the larger fluctuation 
of the ghost propagator in the quenched simulation [65]. 
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Appendix : The gluon propagator, ghost propagator and the running 
coupling in pQCD 


In this appendix, we present the pQCD results of the gluon propagator and the ghost propagator 
and corresponding dressing functions in the MOM scheme, that are used in fitting the lattice 
data. We also present the pQCD definition of the QCD running coupling in the MOM scheme. 
The running coupling of the QCD satisfies the renormalization group equation 


2 


9 


dh 

dq 2 


I3 0 h 2 - /3i/i 3 - foh 4 + ■ ■ ■ 


(A.40) 
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where in general h is scheme and scale dependent. In the high energy region the f3 function of 
the MS scheme is well behaved and expansion in h-^g converges, but in the low energy region, 
expansion in h-^g is not a good converging series. We adopt the MOM scheme and express its 
running coupling in a series of the running coupling of the MS scheme. The method is known as 
the effective charge method|27|. We define as [30] an expansion parameter y^s(q) that satisfies 
asymptotically 

Vs/ms (?) = & log(g 2 /2l|js) - log(/3o3/Mg(g)) ( A - 41 ) 

in terms of the effective coupling in MOM scheme y which is a solution of 

~ = 00 1o &^ 2 / a2 mom) - J 0 lo &(0ov) ( A - 42 ) 

Using the function k(q 2 ,y ) defined as 

Hq 2 ,y) = - + log(/3oj/) - Po log(5 2 / A Ms) (A.43) 

V Po 

the expansion parameter yjjg-(q) is expressed as 

3/ms(?) = J/t 1 + y k (i 2 ,y) + v 2 (tt + Hq 2 , v ) 2 ) 

Po 

+y 3 {^ k (q 2 ,y) + v ) 2 + fc (? 2 > y) 3 ) + ■■■ (A. 44 ) 


In terms of the 3 /^( 5 ), the solution of the renormalization group equation 


f3ol0g h = l + Wo l0giP ° h) 


f h , , 1 Pi 
Jo X X 2 Pox 


Po 


Pox 2 + Pix 3 H-+ P„x n + 2 ' 


(A.45) 


can be expressed as 


Hq) — 2/ms(9)( 4 + VMsil) 2 (0^/0o - (Pi/Po) 2 ) 
+?/ms('?) 3 -{P 3 /P 0 — (Pi/Po) 3 ) H- 


(A.46) 


2 _ 

where /3o = 11- Nf,p\ = 102- —Nf, are scheme independent, and in the MS scheme 

, 2857 5033 Ar 325 „ r2 

^~~2 - W Nf + ~U Nft 

149753 , 1078361 6508 AT 

03 = (—— + 3564C(3) + (- N f - —JV><(3)). 

In the effective charge method, the propagator (—q 2 ) (we use Minkovski metric here) is 

expressed by the scale and scheme invariant propagator D^t(~q 2 ) and a function f(fe') [30[|28] 


DZU-q 2 ) = f(h)D“U-q 2 ). 

Where satisfies the renormalization group equation 

y 2 -j^DZti-q 2 ) = (70/1 + 71/2 2 + 72 h 3 H- )DZt{-q 2 ) 


and 


The general solution of (1A.49I) is 


f(h) = exp j 


dx i{x) 
x P(x) 


(A.47) 

(A.48) 

(A.49) 


f(h) = Ah 70 [1 + (71 - 1oPi)h 


+ — ((7l - 0110) + P 2 I 0 + Pi 70 
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-/?i7i - 027o)h 2 

+(^(7i - 0i7o) 3 + i(7i - 0i7o) 

(72 + 0i 2 7o - 017i - 027o) 

+^(73 - 0i 3 7o + 20J027 o - 037o + 0i 2 7i 
-0a7i “ 0 i72 )h 3 H-] 


where 0- = (i = 1, 2, 3) and 7 ,- = (j = 0,1, 2, 3). 


(A.50) 


A.l Gluon dressing function 

The anomalous dimension 73 of the gluon propagator is 

73 = 73 0 ft + 73ih 2 + 732 ft 3 + 73 3 ft 4 + ' 

where 

13 


2Nf_ 

3 ’ 

59 15 32 

73 1 = 9 y-iV /Y -lV / -, 


73o = Y 


and 


„_,9965 9 

732 _ 27( 288 16^ 3 

+ 18C(3)) + 2N f (- A _ 24C(3)) 

76 2 Nf 44 16 

+—N? + + —A T f 

12 * 3 9 9 * 


, 10596127 1012023 . 8019 

73s = -(-+ ^^<(3) ' -T-C(4) 


768 256 


40905 


C(5) 


+N f ( 


23350603 
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5184 
2 . 43033 

+JV/(- 
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+jV/3( “iii + ? c(3))) 


387649 . . 8955 . . 3355 

- ^5T« 3) + —« 4 > + — 

2017 

-^C(3) - 33C(4)) 


(A.52) 


Using the above gluon held anomalous dimension of four-loop level in MOM scheme [28] 
and the coupling constant h(q) in MS scheme (IA.46I) . which is a function of the parameter y 
defined in the MOM scheme, we derive the gluon propagator Da(—Q 2 ) in the MOM scheme 
as a solution of the renormalization group equation 


9 2 ) = (73 0 ft + 73ift 2 + 73 2 ft 3 + '' ■ )D A (-q 2 ) 

The gluon dressing function Z(q 2 ) = q 2 DA{q 2 ) in Eucledian metric becomes 

3 (l04n 2 - 1974n + 15813) h 


Z = X L k 66~4n [1 _ 


16(33 - 2n ) 2 


+{((l28000n s - 192(53419 + 504<(3))n 4 
+288(1235761 + 5238C(3))n 3 
+108(—56578007 + 772200C(3))n 2 
—324(—153696523 + 6930396C(3))n 
+ 243(—615512003 + 60661656C(3))) h 2 } 


(A.53) 
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/(1536(33 - 2n) 4 ) 

+{(-354549760n 8 - 663552(-79985 + 304C(3))n 7 
+4608(—738019369 + 10620024C(3) + 4370400C(5))n 6 
—3456(—34931893063 + 1008068136C(3) 
+523278720C(5))n 5 

+2592(46615708836^ (3) + 275(-3654988631 
+93932928C(5)))n 4 

—5832(—6037451357147 + 404411943104^(3) 
+223742745600£(5))n 3 

+2916(—100416325711969 + 9138430613136C(3) 
+4824029548800C (5) )n 2 

—4374(—314978703784231 + 37405611077472^(3) 
+18085875033600£(5))n 

+6561(—430343889400537 + 64653527897640C(3) 

+ 27401036762880^(5))) h 3 } 

/(663552(33 - 2n) 6 ) 

+0 (h 4 )} 


where n = Nf. 


A.2 The ghost dressing function 

The anomalous dimension 73 of the ghost propagator is 

73 = 73 0 h + 731 h 2 + 73 2 h 3 + 733 h 4 H - 


where 

73 0 

731 

73 2 = 2 N f (-^ + 12C(3)) + 

9 97 15817 9 

+ 2 iV ^-l08- 9 ^ ) + 27 (W + 32« 3 » 

and 


_ 9 

= 9^ — 31V/ — , 
48 1 12’ 


73 3 = -(- 
40905 


2857419 


512 


19244 ° 7 c ( 3 ) + 8019 c(4) 


+- 


3355 


C(5) + Nf( 
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512 

1239661 

1152 


48857 

48 


64 

C(3) — 


8955 

32 


C(4) 


- f C<3) + f C(4)) 


(A.54) 


Using the ghost field anomalous dimension of four-loop level in MOM scheme and the 
coupling constant h(q) in MOM scheme 1281 . we derive the ghost propagator -Dg(— q 2 ) in the 
MOM scheme as a solution of the renormalization group equation 

Ai 2 ^-C>G(-g 2 ) = (730^ + 731^+732^ H- )D G (-q 2 ) 


(A.55) 
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The ghost dressing function G(q 2 ) = q 2 Dc(q 2 ) in Eucledian metric becomes 

G _1 = A G _1 h“T3^[l + 

( 10 n 3 (40n 2 + 138n - 1611) 97 \ 

^ 9 + 8(33 - 2n) 2 12 J 

+{(512(1439 + 48£(3))n 5 
-3840(13883 + 174C(3))n 4 
—864(—1728454 + 17931C(3))n 3 
+108(—185691691 + 6984516C(3))n 2 
—324(—395301865 + 28831638C(3))n 
+ 19683(—15277259 + 1921964C(3))) h 2 } 

/(1152(33 — 2n) 4 ) 

+{(-16384(174163 + 432C(3))n 8 

— 12288(—30802025 + 637212^(3) + 1002240C(5))n 7 
+4608(—4614333119 + 207142932C(3) 

+266137650C(5))n® 

-3456(-192809757953 + 13588881045C(3) 

+14877266760C(5))n 5 

+5184(—2470563836117 + 240877496568C(3) 

+225684805500C(5))n 4 

— 1944(—79693953595001 + 10028539488648C(3) 
+7956577252800C(5))n 3 

+43740(—26302376345491 + 4087102826048C(3) 
+2675346352272C(5))n 2 

— 13122(—363568314295693 + 67715969212716C(3) 
+34697940156000C(5))n 

+59049(—141629801206331 + 31037533417440C(3) 

+ 11069576361360^(5))) ft 3 } 

/(746496(33 - 2n) 6 ) 

+Q (ft 4 )]. (A.56) 


A.3 The running coupling 

In perturbative QCD, running coupling is derived by the renormalization group equation 


da 

dlogg 


(00 2 . 01 3 . 02 4 

(- CX- H-rtt H- ttCX + 

' 27t 47t 2 647t 3 


03 


128tt 4 


a 5 ) + o(a 6 ) 


(A.57) 


In the MOM scheme, inversion of the 2-loop formula 

+ = Aie “3§fc(^^)- /5l//3 o (A.58) 

47r 

can be done analytically with use of the Lambert W function that satisfies 2 = W[z]e w ^. We 
find 

^ = ^W[(^/2/?i)e (/3 o /2/3l)t ] (A.59) 

i7T +1 

where t = log {g 2 /A 2 ). 

An approximate inversion of the four-loop formula yields the running coupling as a function 
of t = log (g 2 /A 2 ) as follows [36] [39] . 
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%s,pert 


(m) = 


47T 8-7T/3i log(t) 


Pot 


+ 


1 ( 2np 2 

(Pot) 3 V Po 
2np 3 

Po 


(Pot) 4 


+ 


+ 


Po (Pot) 2 
16 nPf 


Po 

I67 rp 3 


(log (t) - log(f) - 1) 


(—2 log 3 (t) + 5 log 2 (t) 


+(4-^)]og(t)-l 


(A.60) 
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